Finite temperature dynamical DMRG and the Drude weight of spin-1/2 chains 
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We propose an easily implemented approach to study time-dependent correlation functions of 
one dimensional systems at finite temperature T using the density matrix renormalization group. 
The entanglement growth inherent to any time-dependent calculation is significantly reduced if the 
auxiliary degrees of freedom which purify the statistical operator are time evolved with the physical 
Hamiltonian but reversed time. We exploit this to investigate the long time behavior of current 
correlation functions of the XXZ spin-1/2 Heisenberg chain. This allows a direct extraction of the 
Drude weight D at intermediate to large T. We find that D is nonzero - and thus transport is 
dissipationless - everywhere in the gapless phase. At low temperatures we establish an upper bound 
to D by comparing with bosonization. 
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It is an intriguing question if a physical system can ex- 
hibit dissipationless transport. In this case the conduc- 
tivity has a singular contribution DS(uj) where D is typi- 
cally referred to as the Drude weight. As a consequence, 
a fraction of an initially excited current will survive to 
infinite time. If the current operator is conserved by the 
Hamiltonian, the corresponding quantum system clearly 
supports dissipationless transport at any temperature T. 
The more subtle question of whether the Drude weight 
can be nonzero when the current operator has no overlap 
with any local conserved quantity has attracted consid- 
erable attention [H4l?i| without being resolved. In one 
spatial dimension it is believed that D ^ is only pos- 
sible at T > for an integrable model where an infinite 
set of conserved local operators might restrict dissipation 
processes. 

This paper uses a finite-temperature time-dependent 
density matrix renormalization group (DMRG) approach 
to calculate the Drude weight for a prototypical inte- 
grable one-dimensional system - the antiferromagnetic 
XXZ Heisenberg chain. The latter describes L — > oo in- 
teracting spin-1/2 degrees of freedom S^' y ' z 



H = jJ2(s:s* +1 + sy n s y 
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or equivalcntly spinless fermions through the Jordan- 
Wigncr transformation. The spectrum is gapless for all 
anisotropics < A < 1. While Bethe ansatz allows eval- 
uation of D at T = via the Kohn formula [lil . [2(J , it 
has not been possible to reliably compute its value for 
nonzero temperatures. The Mazur lower bound 21[ de- 
termined by all local conserved quantities is zero in ab- 
sence of a magnetic field (however, a lower bound was 
recently constructed TfJ 17 1 from a nonlocal operator at 
T = oo; we will explicitly investigate whether or not it 
saturates the Drude weight). Two different Bethe ansatz 



approaches 



sites, and quantum Monte Carlo 0,0 data requires an ill- 
controlled analytic continuation in order to extract D. A 
recent bosonization study combined with DMRG numer- 
ics [HI 111 proves the existence of a diffusive transport 



channel 22j and yields an upper bound on the Drude 
weight; but the timescales reachable within the DMRG 
are yet too small to make a decisive statement about 
whether or not D is nonzero. While most of these works 
conclude in favor of a finite Drude weight in the gapless 
phase, no ultimately accepted quantitative picture (e.g., 
concerning its A and T - dependence and whether or not 
it vanishes for the isotropic chain) has emerged. It is the 
first goal of this work to fill this gap at intermediate and 
high temperatures. 

The density matrix renormalization group p3j was ini- 
tially devised as a powerful tool to calculate ground state 
properties of one-dimensional systems with short-ranged 
interactions 0, dBf and su bseq uently extended to ad- 
dress real-time dynamics [26|-|29(. In princi ple, it can be 
applied straightforwardly to nonzero T [30l - l32j by intro- 
ducing auxiliary degrees of freedom to purify the ther- 
mal statistical operator (33|. In practice, however, the 
increase of entanglement with time has limited previous 
finite-temperature calculations to rather short timescales 
USUI. The second aim of our paper is to propose an 



yield contradictory results. Exact diag- 



onalization [3j-|6| can only treat systems of up to L ~ 20 



easy-to-implement modification to the DMRG algorithm 
of Ref. l3ll : The auxiliaries are by construction inert, but 
they can be exposed to an arbitrary unitary transforma- 
tion without involving any approximation. It turns out 
that the intuitive choice of time-evolving with the phys- 
ical Hamiltonian but reversed time leads to a drastic re- 
duction of the entanglement growth - thus, significantly 
longer timescales can be reached and the Drude weight 
of the XXZ chain can be calculated after all. 

This paper is organized as follows. We first explain our 
modified DMRG method and test it for the exactly solv- 
able case A = 0. It is then used to compute the Drude 
weight for a range of A and T. In a nutshell, we find 
that D is nonzero everywhere from the non-interacting 
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FIG. 1. DMRG calculation for the finite-temperature longi- 
tudinal spin structure factor of a L = 100 site XX Heisen- 
berg chain (A = 0) in comparison with the exact result in 
the thermodynamic limit [34| . We show both the real space 
representation S^ z (t) = {S^^ 2 +n(^)^L/2) an d Fourier trans- 
form S?(t) = E n e lfc ™Sn Z (t). In (a) [in (b) and (c)], the 
DMRG Hilbert space dimension \ [the discarded weight] is 
fixed. Whereas the traditional DMRG approach [3l[ breaks 
down rapidly, larger times can be reached if the auxiliaries 
are evolved with the physical Hamiltonian but reversed time. 



point A = to the isotropic chain A = 1 and decreases 
with increasing temperature and A; values in the gapped 
regime A > 1 are numerically consistent with D = 0. We 
stress that our approach does not involve any approxi- 
mation: The Drude weight can be read off directly from 
the long-time asymptotics of the spin current correlation 
function (l3j ; the system size is large enough (typically 
L ~ 100 — 250) to be in the thermodynamic limit for each 
temperature at hand; and DMRG is numerically exact. 

Finite-temperature DMRG — In order to eventually 
compute correlation functions of the type 5^(t) = 
(Sji(t)S^) = Tr [p SH(t)S^] within the DMRG, one 
needs to purify the thermal grandcanonical density ma- 
trix [3^ by introducing an auxiliary Hilbert space Q: 
pp = Tr q|^)(^|- This is analytically possible only 
at infinite temperature j3 = 1/T = where po = 2~ L . 
However, can be obtained from \^o) by applying 

an imaginary time evolution, \&p) = e~^ H ' 2 \^o) [23| . 
Any correlation function can therefore be exactly re- 
cast as S^(t) = ^p\e iHt S^- im S^p)/y/(^pkp). 
These objects, however, are directly accessible within 
a standard time-dependent DMRG framework 
l28j |. It is most convenient to first express 
in terms of a matrix product state 
Y ja A ai A° 2 ■ ■ ■ A a2L \<7\<J2 ■ ■ ■ 02l) , where 
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36], 



l*o) 
denotes 



the eigenbasis of S* , and Q is spanned by the states with 
even indices. We choose A^ odd — (1 0), A^ odd = (0 — 1), 
Steven = (o l/y/2) T , A±°™ = (1/^/2 0) T in order to ex- 



ploit ^-conservation [23j. After factorizing the imagi- 
nary and real time evolutions operators exp(~ H) using a 
fourth order Trotter decomposition [23[ , they can be suc- 
cessively applied to l^o)- At each time step, two singular 
value decompositions are carried out to update three con- 
secutive matrices A" n (the physical sites are the odd ones 
and thus next-nearest neighbors). The matrix dimension 
X is dynamically increased such that at each time step 
the sum of all squared discarded singular values is kept 
below a threshold value e. The DMRG approximation to 
Snm(t) thus becomes exact in the limit e-)0. 

Efficient DMRG approach to T ^ — Prior finite- 
temperature studies employing the DMRG were in prac- 
tice limited to fairly intermediate timescales due to a 
rapid blow-up of the block Hilbert space dimension \ 



3ll . |34| . While in a ground state calculation the entan- 



glement entropy grows locally around a quench <S^|gs) 
at T > it increases homogeneously even when only the 
thermal density matrix (i.e., the state IvE^)) is exposed 
to a supposedly trivial time evolution. However, the pu- 
rification is not unique - we have the analytic freedom 
to apply any unitary transformation U aux (t) to the in 
principle inert auxiliary sites (physical quantities are de- 
termined by the trace over Q and are thus not affected by 
C4ux)- Thus, it is reasonable to investigate: Is there some 
unitary operator which can undo the above damage? In- 
deed, the intuitive choice U aux (t) = e +lHt , with H being 
the time- reversed physical Hamiltonian applied to the 
auxilliary sites renders the time evolution of \^ p) trivial 
(37j . and the computation of e~ iIrt S^l^ p) is eventually 
only plagued by an entanglement building up around site 
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FIG. 2. Spin current correlation functions of the XXZ Heisen- 
berg chain at intermediate to large temperatures. The long- 
time asymptotics determine the Drude weight D. Inset: 
Bosonization result of Ref. [ll| evaluated for D — 0.05. De- 
spite the large temperature, bosonization agrees qualitatively 
with the DMRG data. 
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FIG. 3. Comprehensive picture of the A-dependence of current correlation function at intermediate (T/J = 0.5) and high 
(T = oo ) temperatures [4[|. The corresponding Drude weight D is shown in the insets and for infinite T compared to the 
lower bound of Ref. Il6l . The left inset to (b) shows data for various system sizes L = 30 . . . 200, illustrating that our DMRG 
calculation can readily reach the thermodynamic limit. 



n [38[. This leads to a significantly slower increase of x, 
and thus longer timescales can be reached. 

To demonstrate the potential of this modification - 
which does not involve any approximation except for the 
usual truncation inherent to each DMRG update - we 
calculate the longitudinal spin structure factor of the XX 
chain and compare with the exact solution obtained by 
mapping to free fermions. The result is shown in Fig. [T] 
Even for times t J ~ 30, x omv grows to moderate x ~ 
250 in order to keep the total (summed) discarded weight 
below e ~ 10~ 6 . This holds for the whole temperature 
range T/J = 0.1 . . . oo and is to be contrasted with the 



DMRG approaches of Refs. |31| and |34| which break down 
much earlier [see Fig. [lja) and note that increasing x 
from 60 to 200 merely allows t J ~ 9 instead of t J ~ 7]. 

In order to assure that the above scenario is generic, 
we have also studied (i) longitudinal and transverse spin 
structure factors of the XX chain in presence of a mag- 
netic field, (ii) the XXZ chain, (iii) the anisotropic spin 
1 Heisenberg chain, and (iv) the quantum Ising model 
and compared with available analytics or DMRG data 
3l[ |34| . The reachable timescales depend on the model 
and correlation function under investigation - but in- 
troducing t/ aux generally reduces the effort of finite- 
temperature DMRG to that of a ground state calcula- 
tion. 

Drude weight at intermediate to large T — We will now 
employ our modified DMRG to extract the Drude weight 
of the XXZ chain in the critical phase < A < 1. The 
spin current is defined from a continuity equation and 



reads j = j Tl 



lL 

2 



J2n( S t S n+i - s n+i s n ), where 



5„ = S% ± iS% (and n denote physical sites only). The 
long-time asymptotics of the current correlation function 
directly yield D [lj]: 



D 



lim 

t— >oo 



2LT 



lim 



2T 



(2) 



A similar strategy was pursued in Ref. [lCj using the trans- 
fer matrix DMRG. We quantitatively reproduce those re- 
sults but are able to reach time regimes about a factor 
of 2 — 3 larger within one day of computer time (39j. At 
intermediate temperatures T/J > 0.5, this allows us to 
access a scale where (j(t)j) is saturated, which enables a 
quantitative estimate of the Drude weight. This is illus- 
trated in Fig.[2]for fixed A = 0.6 and three temperatures 
T/J = 0.5,2,oo. The current correlators become fiat 
around U ~ 20, and we obtain D ~ 0.06J/2T for all 
T/J > 0.5. Note that at T/J = 0.5, the Drude weight is 
reduced by a factor of 3 compared to its zero-temperature 
Bethe ansatz value [HI . We emphasize that for each tem- 
perature we chose a system size large enough to deter- 
mine D in the thermodynamic limit [this is demonstrated 
explicitly in the inset to Fig.[3Jb)]. 

It is instructive to compare our data with the bosoniza- 
tion calculation of Ref. [ll| which yields the current corre- 
lator up to the asymptotic value D. With our prediction 
for D, the bosonization curve at T/J = 0.5 agrees quali- 
tatively with the DMRG curve despite the large tempera- 
ture (see the inset to Fig. [2]) . This demonstrates that it is 
indeed reasonable to expect (j(t)j) to saturate at a scale 
t J ~ 20. At larger temperatures, the correlators evolve 
non- monotonously at intermediate t [Io| . The magnitude 
of these oscillations, however, decreases as T is lowered; 
they eventually die out around T/J = 0.5. 

We now investigate the A-dependence of the Drude 
weight for intermediate (T/J = 0.5) and high (T = oo) 
temperatures. As illustrated in Fig. D monotonously 
decreases from its trivial value at A = (where j com- 
mutes with the Hamiltonian and thus (j(t)j) = const.) 
down to a finite value for the isotropic chain. The latter 
is consistent with most previous works which event ually 
concluded in favor of D(A = 1) > @, |, i, 1, (Tj. 
The qualitative behavior -D(Ai) < D(A 2 ) for Ai > A 2 
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FIG. 4. Spin current correlation functions at low tempera- 
tures. Due to the presence of a diffusive transport channel [l(| 
with a small decay rate 7 (manifesting as a linear decrease of 
(j(t)j) for times 1/T Ct< I/7), the Drude weight can only 
be estimated indirectly by comparing with the bosonization 
prediction [HI]. For T/J = 0.2 and A/J = 0.6, this suggests 
D ~ 0.07 as an upper bound. 



agrees with one [l[ of the Bethe ansatz [l], Q calcula- 
tions [note, however, that Ref. [l] finds D(A = 1) = 
for any nonzero T]. At infinite temperature, our data is 
consistent with a recently-established lower Mazur bound 
[lij . The latter has a fractal A-dependence and seems 
to saturate our Drude weight for commensurate values 
A = cos(ttI / m) , I , m € Z + [such as A = 0.5; see the left 
inset to Fig.^b)]. 

Low temperatures — At low T, bosonization [ToL HH 
conjectures a diffusive transport contribution: (j(t)j) 
falls off exponentially for times t > 1/T and saturates 
at a scale set by the inverse decay rate 7 (bosonization 
does not yield the saturation value ~ D). This picture 
is strongly supported by DMRG [Io[ as well as Quan- 



tum Monte Carlo 221 numerics. Our modified DMRG 



reconfirms the existence of the diffusive channel, and we 
extract 7 (by fitting the exponential in its linear regime 
t -C I/7) in nice agreement with the bosonization predic- 
tion (see Fig. 2] keep in mind that the DMRG calculation 
of Ref. [13 reached tJ ~ 7). Unfortunately, 7 decreases 
with temperature and for T/J < 0.4 we can no longer 
access the scale t » 1/7 where the current correlator 
saturates. 

We will now discuss an indirect way to determine D 
which was proposed in Ref. 111! The only free parameter 
in this bosonization calculation for (j(t)j) is the Drude 
weight, which one can therefore try to extract by fitting 
the DMRG data at intermediate times. This is exem- 
plified in Fig. H for A = 0.6 and T/J = 0.2. Our nu- 
merics suggests D ~ 0.07 as an upper bound (and seem- 
ingly contradicts both Bethe ansatz predictions D = 0.15 
Q and D = 0.105 This in turn indicates that 

the Drude weight only slightly deviates from its value 
D ~ 0.06 at T/J = 0.5 when temperature is lowered 



to T/J = 0.2, and thus D supposedly varies strongly 
below T/J < 0.2 in order to eventually approach its 
zero-temperature value D = 0.15 [p|. We emphasize, 
however, that these arguments are less stringent than 
the direct observation of a flattening (j(t)j): T/J = 0.2 
might be too high for the low-energy bosonization ap- 
proach to be ultimately accurate - but at lower T, we 
cannot access the time regime 1/T < t < 1/7 necessary 
for a meaningful comparison. 

Gapped phase — The gapped phase of the Heisenberg 
chain has attracted less attention than the critical one: 
At T = 0, the Drude weight vanishes for all A > 1 [l9j], 
so one might reasonably expect the same to hold at finite 
temperature [f| [l(| EE EH- Indeed, our data for T/J = 
0.5 and A = 1.2, 1.4, 1.6 is consistent with a zero Drude 
weight [sec Fig.[3{a)]. However, the effort to reach large 
times increases significantly with A, and we thus leave a 
detailed investigation of D in the close vicinity of A > 1 
subject for future work. 

Conclusion and Outlook — We have introduced a mod- 
ification to finite-temperature DMRG which can be eas- 
ily incorporated within existing implementations of the 
method. It allows to compute correlation functions up 
to significantly larger times. We have exploited this to 
quantitatively extract the Drude weight D of the antifer- 
romagnetic spin 1/2 XXZ Heisenberg chain at interme- 
diate to large temperatures. D decreases monotonously 
when the asymmetry parameter A is varied from A = 
(XX chain) to A = 1 (isotropic chain). Our data strongly 
suggests that D is finite at A = 1. For small T, we recon- 
firmed the existence of a diffusive transport contribution 
in the critical phase. Our access to the Drude weight, 
however, is still only qualitative at low temperatures, 
and quantitatively bridging the gap between zero tem- 
perature (where one can carry out Bethe ansatz calcula- 
tions) and our data at intermediate and large T remains 
as the major open challenge. The DMRG algorithm pre- 
sented may enable various applications - investigation of 
spin and current correlation functions of non-intcgrablc 
models being an obvious one. 
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